function fig3_15

% exact solution of ep*y"+(x-a)*(x-b)*(4y'-1)=0

global a b ep

a=0.249
%a=0.252

nx=4000;
y(1)=-2;
b=0.75
ep=0.0002;

x=linspace(0,1,nx);
for i=2:nx
    dy=int2(@f,x(i-1),x(i));
    if i==2
        y(i)=dy;
    else
        y(i)=y(i-1)+dy;
    end
end
A=(5-0.25)/y(nx);
for i=2:nx
    y(i)=-2+0.25*x(i)+A*y(i);
end

clf
% get(gcf)
set(gcf,'Position', [25 1115 658 230])
co = [0 0 1;
    0 0.5 0;
    1 0 0;
    0 0.75 0.75;
    0.75 0 0.75;
    0.75 0.75 0;
    0.25 0.25 0.25];
set(groot,'defaultAxesColorOrder',co)
subaxis(1,1,1,1,'MT',0.01,'MB',0.18,'MR',-0.02,'ML',0.05,'P',0.04)

plot(x,y,'b','LineWidth',1.3)
box on
grid on
axis([-0.01 1 -2 3.1])

set(gca,'ytick',[-2  0  2])
set(gca,'YTickLabel',{'$-2$';'0';'2'})

ax = gca;
ax.XAxis.TickLabelInterpreter = 'latex';
ax.YAxis.TickLabelInterpreter = 'latex';
set(gca,'FontSize',13)

say=['$a = \,\,$',num2str(a)];
text(0.4,-0.5,say,'FontSize',14,'FontWeight','bold','Interpreter','latex')

xlabel('$x$\,-\,axis','Interpreter','latex','FontSize',14)
ylabel('Solution','Interpreter','latex','FontSize',14)

% if a==0.252
%     exportgraphics(gcf,'/Users/mark/Desktop/a1.pdf')
% elseif a==0.249
%     exportgraphics(gcf,'/Users/mark/Desktop/a2.pdf')
% end

function g=f(x)
global a b ep
g=exp(-4*(x^3/3-0.5*(a+b)*x^2+a*b*x)/ep);


